subroutine calibrate
use params
use mini
implicit none
include 'mpif.h'

double precision :: temp

double precision, parameter :: limit_2004=0.537461728, pop_debt_2004=0.301624055, access_young_2004=.5792
double precision, parameter :: acc_2004=0.716222961, def_unmodeled_data=0.448051948, SS_2004=0.295732318
double precision, parameter :: credit_2004=0.079282605, def_rate_2004=0.003798383
double precision, parameter :: limit_2016=0.381465888, def_rate_2019=0.0014
double precision, parameter :: tol_cal=1e-3

double precision :: avg_SS
double precision, parameter :: credit_drop_data=0.023153, limit_relative_drop_data=.2294
double precision :: credit_drop_model, limit_relative_drop_model
double precision :: d_cal
integer, parameter :: itmax_calibrate=1

double precision, parameter :: pow=.75

integer, parameter :: n_amoeba=4
integer :: iter_amoeba
double precision :: tol_amoeba=1e-4
double precision :: x_amoeba(n_amoeba), y_amoeba, x_optimal(n_amoeba), y_optimal
integer :: find, find_amoeba
double precision :: x_amoeba_recmpi(n_amoeba*nproc), y_amoeba_recmpi(nproc)

! if ( rank .eq. 0 ) then

 allocate (rn(itmax_calibrate,nproc,n_amoeba+1,n_amoeba))
 
    if ( rank .eq. 0 ) then
    
        do iter_calibrate=1,itmax_calibrate
        do iproc=1,nproc
        do iter=1,n_amoeba+1
        
            call random_number(rn(iter_calibrate,iproc,iter,1))

            call random_number(rn(iter_calibrate,iproc,iter,2))

            call random_number(rn(iter_calibrate,iproc,iter,3))

            call random_number(rn(iter_calibrate,iproc,iter,4))
           
        end do
        end do        
        end do
        
    end if

    do iter_calibrate=1,itmax_calibrate
    do iproc=1,nproc
    do iter=1,n_amoeba+1
   
        call MPI_BCAST (rn(iter_calibrate,iproc,iter,:),4,mpi_double_precision,0,mpi_comm_world,ierr)
    
    end do
    end do
    end do
    
var_z=.063
var_eta=.024
var_nu=.244
rho_eta=.977
    
    if ( rank .eq. 0 ) then
    print*,''
    print*,'original'
    print*,'transitory_variance'
    write(6,'(300F14.10)')var_z  
    print*,'persistent_variance'
    write(6,'(300F14.10)')var_eta
    print*,'permanent_variance'
    write(6,'(300F14.10)')var_nu
    print*,'persistence'
    write(6,'(300F14.10)')rho_eta
    print*,''
    end if
    
     open (unit=1,file='input_parameters.txt')
     read(1,*)x_optimal
     close(1)

y_optimal=func_sum(x_optimal)
find=1
        
    do iter_calibrate=1,itmax_calibrate
    
        call amoeba (x_amoeba,y_amoeba,n_amoeba,tol_amoeba,iter_amoeba,func_sum)

        call MPI_ALLGATHER (x_amoeba,n_amoeba,mpi_double_precision,x_amoeba_recmpi,n_amoeba,mpi_double_precision,mpi_comm_world,ierr)
        
        call MPI_ALLGATHER (y_amoeba,1,mpi_double_precision,y_amoeba_recmpi,1,mpi_double_precision,mpi_comm_world,ierr)

    find_amoeba=0
    
        do iter=1,nproc
        
            if ( find_amoeba .eq. 0 ) then
            
            find_amoeba=1
            y_amoeba=y_amoeba_recmpi(iter)
            x_amoeba(1)=x_amoeba_recmpi((iter-1)*4+1)
            x_amoeba(2)=x_amoeba_recmpi((iter-1)*4+2)
            x_amoeba(3)=x_amoeba_recmpi((iter-1)*4+3)
            x_amoeba(4)=x_amoeba_recmpi((iter-1)*4+4)

            else if ( find_amoeba .eq. 1 ) then
            
                if ( y_amoeba_recmpi(iter) .le. y_amoeba ) then
                y_amoeba=y_amoeba_recmpi(iter)
                x_amoeba(1)=x_amoeba_recmpi((iter-1)*4+1)
                x_amoeba(2)=x_amoeba_recmpi((iter-1)*4+2)
                x_amoeba(3)=x_amoeba_recmpi((iter-1)*4+3)
                x_amoeba(4)=x_amoeba_recmpi((iter-1)*4+4)
                end if
            
            end if
            
        end do
            
        if ( find .eq. 0 ) then
        
        find=1
        x_optimal=x_amoeba
        y_optimal=y_amoeba
        
        else if ( find .eq. 1 ) then
        
            if ( y_amoeba .le. y_optimal ) then
            x_optimal=x_amoeba
            y_optimal=y_amoeba
            end if
            
        end if        
        
        if ( rank .eq. 0 ) then
        
        print*,''
        print*,'iter_calibrate',iter_calibrate
        print*,'var_z var_eta var_nu rho_eta y_optimal'
        write(6,'(300F14.10)')x_optimal,y_optimal

        var_z=x_optimal(1)
        var_eta=x_optimal(2)
        var_nu=x_optimal(3)
        rho_eta=x_optimal(4)

        var_z=x_optimal(1)
        var_eta=x_optimal(2)
        var_nu=x_optimal(3)
        rho_eta=x_optimal(4)
    
            call gridearnings

        print*,'std_log_earnings_change_1 data model'
        write(6,'(300F14.10)')std_log_earnings_change_1_2019,std_log_earnings_change_1
        print*,'std_log_earnings_change_5 data model'
        write(6,'(300F14.10)')std_log_earnings_change_5_2019,std_log_earnings_change_5
        print*,'std_log_earnings_permanent data model'
        write(6,'(300F14.10)')std_log_earnings_perm_2019,std_log_earnings_perm
        print*,'correlation_lag_1 data model'
        write(6,'(300F14.10)')correlation_lag_1_2019,correlation_lag_t(1)
        print*,'correlation_lag_2 data model'
        write(6,'(300F14.10)')correlation_lag_2_2019,correlation_lag_t(2)
        print*,'correlation_lag_3 data model'
        write(6,'(300F14.10)')correlation_lag_3_2019,correlation_lag_t(3)
        print*,'correlation_lag_4 data model'
        write(6,'(300F14.10)')correlation_lag_4_2019,correlation_lag_t(4)
        print*,'correlation_lag_5 data model'
        write(6,'(300F14.10)')correlation_lag_5_2019,correlation_lag_t(5)
    
            if ( isnan(var_z) .or. isnan(var_eta) .or. isnan(var_nu) .or. isnan(rho_eta) ) then
            print*,'isnan'
                call exit
            end if
    
!            open (unit=1,file='input_parameters.txt')
!            write(1,*)var_z,var_eta,var_nu,rho_eta
!            close(1)
        
        end if
        
    end do    

 deallocate (rn)
! end if

end subroutine calibrate
